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Rocks were subjected to the deformation test under five 
different confining stresses (0, 5, 10, 15, and 20 MPa) using 
the Hoek cell to determine changes in the elastic properties 
of the rocks under confining stress, and the results were 
evaluated based on density, porosity, Schmidt hardness, and 
compressive strength. A total of nine different rocks, two 
granites, two andesites, two limestones, one tuff, one diorite, 
and one marble, were used. When the confining stress was 
increased from 0 MPa to 5 MPa and from 5 MPa to 10 MPa, 
elasticity increased by approximately 20%. When the 
confining stress was increased from 10 MPa to 20 MPa, the 
increase was 7% in comparison with the previous value. 
Then, to formulate the modulus of elasticity for rocks under 
the triaxial stress conditions, a new and intelligent approach 
to gene application, gene expression programming was 
utilized. The success of the model was thoroughly assessed 
based on measurable criteria such as the root mean square 
error, mean absolute percentage error, and coefficient of 
determination. Furthermore, the success of the model was 
comprehensively assessed based on the model testing, and 
0.88 and 0.81 R° values were obtained for training and 
validation, respectively. The performance of the gene 
expression programming-based formulation was compared 
with the formulae previously proposed in the literature. The 
gene expression method exhibited the best performance, and 
it was identified to calculate the modulus of elasticity under 
triaxial stress conditions more effectively. 
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1. Introduction 


The characterization of the engineering properties of rocks plays a vital role in planning the 
engineering structure in terms of stability and cost analysis. The modulus of elasticity (E) 
represents an essential parameter for understanding the deformation behavior of rocks and is 
among the most significant mechanical properties of rocks concerning their areas of usage. The 
said parameter is crucial in tunnel projects, slope consistency, pillar configuration, rock 
destruction and drilling, dams, and numerous other civil and mining applications. It has been 
widely utilized to analyze structural deformations, shrinkage, creep, crack control, etc. 


The modulus of elasticity is used for the purpose of determining the deformation properties of 
wood, metal, concrete, rocks, and many other materials. The mentioned parameter represents the 
ratio of stress to corresponding strain in material under tension or compression. It may be 
obtained from unconfined compressive strength (G;) tests and the incline of tensile test stress- 
strain diagrams for rocks. The deformation test performed to determine the modulus of elasticity 
is difficult, time-consuming, and requires particular specimens and equipment. Many authors 
have tried to identify this parameter using more easily identified material properties due to the 
difficulties in determining the modulus of elasticity using regression analysis [1-4], genetic 
programming applications, e.g. artificial neural networks [4-8], fuzzy logic [9-11], and genetic 
algorithms [12-14]. 


The elasticity test is performed only under uniaxial stresses. Nevertheless, in practice, 
underground rocks are also affected by confining stress (G3), assuming that G2= G340 and the 
modulus of elasticity (E) may be insufficient to reflect the actual conditions of rocks. In real 
situations, rocks are exposed to confining stress, and the effects of their elasticity characteristics 
under confining stress vary. In previous studies, it was stated that the strength and the modulus of 
elasticity of rocks would increase with the increase in G3. However, there are not many studies in 
the literature on the rate of this increase and estimating the respective values. Liang et al. [15] 
and Ma et al. [16] investigated the stress-strain relationships for rock salt and anhydrite. On the 
contrary, Li et al. [17] reported that the modulus of elasticity was independent of strain rate. 
Hoek and Brown [18], Arora [19], and Asef and Reddish [20] developed a different formula for 
predicting the modulus of elasticity with confining pressure. 


This study's goal was to predict E at a triaxial stress state using a novel approach, gene 
expression programming (GEP). The complicated relationship between the parameters affecting 
the E value may be easily modeled using the GEP approach, unlike classical methods. Two 
hundred fifty-six data sets were collected from experimental studies using nine different rocks 
(two granites, two andesites, two limestones, one tuff, one diorite, and one marble). E was used 
as an output, and five parameters, including G3 (MPa), G; (MPa), porosity (g/em*) (8), ultrasonic 
sound velocity (km/s), and Schmidt rebound hardness (SH), were used as the input parameters. 
The model was trained with 192 data sets of the experimental results, following which 64 data 
sets were utilized for testing. The result of the obtained GEP-based formulation was tested and 
compared to the formulations proposed by Arora [19] and Asef and Reddish [20] using 
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experimental data. It was observed that the findings acquired using the model corresponded to 
the experimental data at a high rate. 


2. Material and experimental methods 


2.1. Rock characteristics 


The physical characteristics of the sampled rocks were identified by @ and ultrasonic P-wave 
velocity. The samples were cored from large blocks in the laboratory, and their testing was 
carried out following the procedures described in ISRM [21] . Density (9) and 6 were found by 
employing the water saturation and caliper techniques on the NX-sized samples, using five 
samples every time. The high-frequency ultrasonic pulse technique described by ASTM [22] was 
employed to measure the velocities of the P (compression) waves using the 1 MHz nominal 
frequency. The wave velocity through the sample was computed using the travel time from the 
transmitter to a receiver. A minimum of five measurements was made from every rock type, and 
they were averaged. 


The measured mechanical properties were UCS and rebound tests, which were determined in 
accordance with the standards by ISRM [21]. A minimum of five core samples from every rock 
were exposed to strength tests, which were conducted using a fully automatic, servo-controlled, 
instrumented, and computer-controlled press machine. Rebound tests were conducted on the 
block samples using the Digi-Schmidt 2000 concrete test hammer. At least twenty measurements 
at various points of every sample were performed for the rebound hardness of rocks. Table 1 
summarizes the physical and mechanical characteristics of the rocks examined. This table shows 
that the standard deviation values were high because of working with natural materials. 


Table 1 
Basic physical and mechanical characteristics of rocks. 


Min. Max. Mean Std. Dev. 


6 (%) 0.36 18.23 4.36 6.16 
p(gricm’) 1.85 2.86 2.57 0.33 
Vp(Km/s) 2.17 5.82 4.28 1.38 
SH 43.00 70.60 58.18 8.05 


G, (MPa) 28.00 167.97 81.14 34.95 
E (MPa) 1278.20 28260.20 10854.07 5513.98 


2.2. Elasticity at a triaxial stress state 


Triaxial compression tests with different G3 were carried out according to ISRM [21] using a 
hydraulic servo-control test system comprised of three major sub-systems, including a servo- 
controlled press machine, automatic axially loading unit, instrumentation, and a personal 
computer. G3 was achieved through hydrostatic pressure using an automatic pressure system for 
confining stress in a Hoek triaxial cell. A series of G3 at 5, 10, and 15 MPa and at 20 MPa for 
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only two rocks (one limestone and one granite) were applied as confining pressure. The 
maximum confining stress was limited to 20 MPa, which was equivalent to a load of 
approximately 800 m depth representing the limit of the depth of most engineering works. Axial 
loading and G3 were simultaneously increased to the pre-established values, following which G3 
was fixed, and axial loading continued to be increased at a constant speed up to the pre- 
determined stresses. The axial strains of the test samples were measured by the linear variable 
differential transformer (LVDT) using the software prepared for controlling the press machine. A 
PC continuously showed and performed the recording of the axial load and the axial 
displacement during the test. The tangential E was derived as a tangent to the stress-strain curve 
(E=AG/Ae) at 50% of the maximum stress (G;) derived (Fig. 1). The E values for different G3 
values of the nine various rocks used in the study are presented in Fig. 2. It is observed that when 
G3 was increased from 0 to 5 and from 5 to 10, the increase in the modulus of elasticity was 
20%, and when G3 was increased from 10 to 15, the mean modulus of elasticity increased by 7% 
on average in the used rocks. 
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Fig. 1. Determination of the modulus of elasticity. 
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Fig. 2. The correlation between E and G3. 


The effects of G; and G3 on the E values are shown in Fig. 3. In particular, when G3 was applied 
as 5 and 10 MPa, it was observed that the modulus of elasticity was about 20% higher compared 
to the previous value. When G3 was applied as 20 MPa, the modulus of elasticity continued to 
increase. However, this rate was only 7% higher than the value at 15 MPa. The most important 
reason for this was the fact that the discontinuities in the studied rocks were continuous to trap 
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up to 10 MPa. After 10 MPa, the closing procedure was continued, but the effect on deformation 
was reduced. In case of the low confining stress, the rocks with low G; and high @ regarding the 
increase in the modulus of elasticity were relatively low. When the confining stress increased, the 
increase in the modulus of elasticity was higher. However, in the rocks with high G; and low 6, 
the situation was the opposite (Fig. 4). Likewise, the highest rate of increase in the modulus of 
elasticity was observed in rocks with the highest SH as the confining stress increased (Fig. 5). 
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Fig. 3. The correlation between E and G;, Gs. 
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Fig. 5. The correlation between E and SH, G3. 
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3. Gene expression programming (GEP) 


Although classical modeling methods such as regression analysis have many advantages, these 
methods need a correlation between independent and dependent variables [23]. The significant 
weakness of traditional modeling methods is their being not equipped for providing useful 
prediction conditions for complicated cases. The soft computing approach surpasses the 
problems and disadvantages of regression analysis and yields incredible and accurate results [9]. 
To overcome this limitation, new methodologies have been proposed. The researchers have 
produced rearranged prediction conditions without accepting an earlier type of the current 
relationship. Gene expression programming (GEP) is another amazing, transformative Al-based 
technique created by Ferreira [24]. It is a groundbreaking developmental calculation that 
combines both the basic direct chromosomes of a settled length like the ones utilized in 
hereditary calculations and the ramified structures of various sizes and shapes like the parse trees 
of hereditary programming. Its assessment arrangement of any information resembles that of the 
biological assessment and is a PC program that is encoded in direct chromosomes of a settled 
length. In the mentioned methodology, a scientific capacity characterized as a chromosome with 
multiple qualities is produced by utilizing the information exhibited to it. 


To create a GEP demonstration, five segments, including the function set, terminal set, control 
parameters, fitness function, and stop condition, are needed. After encoding the issue for the 
competitor arrangement and determining the wellness work, the calculation randomly makes an 
underlying populace of reasonable people (chromosomes) and afterward, changes over every 
chromosome into an articulation tree related to a scientific articulation. A short time later, the 
anticipated target is contrasted, and the real one and the wellness score for every chromosome 
are resolved. If it is adequately high, the calculation stops. Nevertheless, a part of the 
chromosomes are chosen by utilizing roulette wheel testing, and after this process, they are 
transformed to acquire the new ages. This closed cycle proceeds until the wanted wellness score 
is achieved, following which the decoding of the chromosomes is performed for the best 
arrangement of the issue [25,26]. 


GEP has two languages, such as the language of qualities and the language of ETs. In GEP, on 
account of the basic decisions about the structure of ETs and their associations, it is reasonable to 
promptly induce the phenotype presented as the arrangement of quality and vice versa. This 
comprehensible bilingual documentation is named the Karva language [24]. For example, a 
numerical articulation [a * (a/b)| + [‘/a/b] may be expressed by a two-quality chromosome or 
an ET, as seen in Fig. 6. The said figure indicates how two qualities are encoded as a straight 
string and how they are communicated as an ET. The reader is referred to the study conducted by 
Ferreira [24,27] to obtain more detailed information. 


3.1. Overview of gene expression programming 


Due to the failed linear, nonlinear either single or multiple regression in terms of providing the 
suitable outcomes, in this study, the principal point of improvement of GEP-based definitions 
was to produce the numerical capacities for the prediction of E at 5, 10, 15, and 20 MPa G3 at a 
triaxial stress state. While training and testing the GEP model, G3, density (p), 8, SH, Vp, and G; 
were used as the input variables, whereas E was used as the output variable. 
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Fig. 6. The case of the GEP expression tree. 


The GEP model was created with 256 data sets of the rock samples acquired from this 
experimental research. The numbers of the experimental data sets utilized for initial practices and 
testing/validation in the mentioned model were 191 and 63, respectively. Table 2 summarizes the 
parameters utilized in the development of the GEP model. 


Many GEP models were tested to find the best model using the GeneXproTools software. Four 
basic arithmetic operations (+, —, *, /) and several primary mathematical functions (Exp, Ln, x. 
3Rt, Atan, Min, Max, Avg2, Tanh, Not, Inv) were used, and the RMSE limitation was determined 
as the fitness function. 


In the GEP analysis, the setting parameter values significantly affected the output model's fitness. 
They contained the chromosome number, the number of genes, the head size of the gene, and the 
ratio of genetic operators. To select the chromosomal tree, in other words, the head's length and 
the gene number, the GEP models primarily utilized a single gene and two lengths of heads and 
increased the number of genes and heads, consecutively, in every run, while the testing and 
training performance of every model was monitored. Following a number of tests, for the 
models, the numbers of genes 3 and eight heads were revealed to yield the best outcomes. The 
linking of the sub-ETs (genes) was performed as a result of subtraction. 


Table 2 

GEP parameters utilized for the models created. 
Definition of the parameter GEP Model 
Fitness function RMSE 
Number of chromosomes 30 

Gene number 3 

Mutation rate 0.00138 
Inversion rate 0.00546 
One-point recombination rate 0.00277 
Two-point recombination rate 0.00277 
Gene recombination rate 0.00277 
Gene transposition rate 0.00277 
Literals 17 
Arithmetic operators +-x/ 
Mathematical function Exp, Ln, xe. 3Rt, Atan, Min, Max, Avg2, Tanh, Not, Inv 
Head size 8 

Tail Size 9 

Gene Size 26 


Linking function Subtraction 
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Finally, an association of all genetic operators, including transposition, mutation, and crossover, 
was utilized as the set of genetic operators. The chromosome 30 was determined as the best 
generation of individuals in estimating E. The accurate formulation after replacing the constants 
based on the constructed models for E is given by the following formula: 

_ 10"? : 2 = 
pe {68.564 (se “(0% fe) — {min(2094.579¥72; ((Vp + p) * (os * SH)} 


{a (Vp *04 


) « (p* sH)| * ATan(—17.396)} (1) 


Fig. 7 presents the expression tree of the formulation for the model, where dO, dl, d2, and d3 
denote G3, p, 8, SH, Vp, and G;. Table 3 contains the list of the constants utilized in the 
formulation. 


Table 3 
Coefficients used in the GEP model. 
Constant GeneI GenelII Gene III 


CO 8.45 4.09 -2.93 
Cl 6.19 8.65 6.02 
C2 -4.61  -4.79 1.65 
C3 5.49 -204.14 7.46 
C4 2.65 8.53 -2.04 
CS -0.39 — -2.55 -9.48 
C6 -6.03 = -8.69 -17.39 
C7 6.19 -5.13 1.67 
C8 9.13 0.27 9.00 
C9 -9.63 9.27 -5.84 


4. Results and discussion 


The current section investigates the findings of the analysis of the created model and a 
quantitative evaluation of the model's predictive abilities. In assessing the model, defining the 
model's performance and the criteria for estimation accuracy is essential. To produce the GEP 
models and demonstrate the prediction ability of the GEP, the database consisting of the test 
information was separated into two sets, training and testing. From the 256 data sets, for the 
training and testing processes, 192 and 64 experimental data sets, respectively, were selected in a 
random way. 


To compare the performance of the model with the previous models in the literature, three 
previously utilized models were determined. These models are the models of Hoek and Brown 
[18], Arora[19], and Asef and Reddish [20]. 


Hoek and Brown developed a highly accepted formula that uses the modulus of elasticity, strain, 
Poisson's ratio, and confining pressure. This is shown by equation (2). However, it is both very 
expensive and technically difficult to determine the Poisson's ratio used in this equation under 
triaxial stress conditions. 
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Fig. 7. Expression tree of the GEP approach model. 


Ee, = 0, —9(0) + 03) (2) 


Where % is the Poisson’s ratio. 
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Arora suggested a model using three different rock types (plaster of Paris, Jamrani sandstone, 
and Agra sandstone) upon the estimation of the modulus of elasticity with confining pressure. 
This empirical equation is given below: 


E(a3=0) O71 
+ =1- —0.1— 
E(¢3=02%0) ne ( a) (3) 
Asef and Reddish developed the model with its own test data and proposed a new model. It is 
stated that this model illustrates a significant improvement by obtaining the coefficient of 
determination (R7:0.69) with the original data of Arora and its own experimental data using 
arkose sandstone. This empirical equation is presented below: 


Eoszaxe0) _ 200(58) +2 (4) 
E(¢3=0) (23)+p 


Where b=15+60e 2 891 E(g3=0) is the modulus of elasticity at the unconfined compressive 
strength (GPa), E(g3= G2¥0) 1s the modulus of elasticity at a given confining stress (GPa), G3 and 
G; are the confining stress (MPa) and unconfined compressive strength (MPa), respectively. 


The results of the created model utilizing GEP and models proposed by Arora and Ased and 
Reddish were assessed based on the root mean square error (RMSE) values, mean absolute 
percentage error (MAPE), and the coefficient of determination (R’). The above-mentioned 
criteria are presented in equations 5, 6, and 7, respectively. 


_ 1 [SE Iti-oil 
MAPE =~ eae x 100| (5) 


RMSE = ie n (t; —0;)? (6) 


Ree | (n YE troi-(E HEE 01)” | 


7 
(n dR, t?-(2%, ti) )(ndM, 0? -(SL, oi) ) 


Where t denotes the experimental value, o refers to the predicted value, and n refers to the total 
number of data. 


RMSE indicates the sample standard deviation between the observed and predicted values, while 
MAPE represents a measure of the prediction accuracy of the forecasting method in statistics, for 
instance, in the estimation of trends, and it generally presents accuracy as a percentage. R? was 
used as the criterion for evaluating the agreement between the experimental and predicted values. 
The statistical performances of the model established by GEP and the models proposed by Arora 
and Asef and Reddish are summarized in Table 4. As seen in Table 4, for both the training and 
validation sets, the GEP model had low MAPE, RMSE and high R’ values. The R? values related 
to the experimental and predicted data from the GEP model for training and testing were 0.88 
and 0.81, respectively, implying that the model was capable of explaining 88% and 0.81% of the 
variation in the data. 
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Table 4 


Performance statistics of the model. 


GEP Model Arora’s Model Asef and Reddish's Model 

Training Validation Training Validation Training Validation 

MAPE 0.069 0.264 0.173 0.465 1.198 3.80 
RMSE = 1904.84 2688.59 3991.40 4085.46 27206.39 30446.23 
R° 0.88 0.81 0.56 0.60 0,62 0.62 
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The graphical comparison of the E values predicted by the GEP-based models and models 
proposed by Arora and Asef and Reddish in the training and validation phases to their 
experimental values is presented in Figs. 8 and 9, respectively. As is observed from the figures in 
question, the experimental values were close to the estimated values for the GEP-based models. 
They demonstrated the successful performance of GEP models in predicting E at a triaxial stress 
state using 9, p, Vp, SH, and G;. The models proposed by Arora and Asef and Reddish exhibited 
a similar performance to each other in the current experimental data. However, the performances 
of these models are lower than that of the GEP-based model. 


Fig. 8. 


Fig. 9. Comparison of the measured and predicted E results of the validation set. 


160000 


oGEP R? = 0,8828 
140000 ~SArora R*.=.0,5592 A 
Asef and Reddish R* = 0,6205 
120000 A 
© 
= 100000 x 
pty 
> 80000 
£ 
° 
5 60000 
2 
ao 40000 
20000 
0 
0 5000 10000 15000 20000 25000 30000 


Experimental E (MPa) 


Comparison of the measured and predicted E results of the training set. 


140000 -eGEP— Re =0,8105 
oArora. RR? = 0.5995 
120000 aAsefand Reddish R?=0,6206 


100000 
oa 


fe) 
[o) 
fo) 
[o) 
[o) 


Predicted E (M 
ro) 
oO 
oO 
ro) 
ro) 


0 5000 10000 15000 20000 25000 30000 
Experimental E (MPa) 


58 U. Atici/ Journal of Soft Computing in Civil Engineering 4-4 (2020) 47-60 


5. Conclusion 


In this study, the triaxial compression test for nine different rocks (two granites, two andesites, 
two limestones, one tuff, one diorite, and one marble) was carried out at four confining stress 
values (5, 10, 15, and 20 Mpa) to determine the behavior of rocks at a triaxial stress state. The 
GEP-based formulation is displayed for predicting the modulus of elasticity (E) from some mass 
characteristics of rocks by utilizing the GEP system, which is a novel methodology for this kind 
of issues. In the GEP model, confining stress, density, porosity, Schmidt hardness, ultrasonic p- 
wave velocity, and uniaxial compressive strength were utilized as the input variables, whereas E 
was used as the output variable. For this purpose, 256 data sets were used. The comparison of the 
performance of the GEP-based model with the formulae previously suggested in the literature by 
Arora [19] and Asef and Reddish [20] was made using the root mean square error (RMSE) 
values, mean absolute percentage error (MAPE), and the coefficient of determination (R’). 


From the present study, the following was concluded: 


e With the increase in the support pressure, there was a considerable increase in the modulus 
of elasticity. However, when the support pressure was further increased, a reduction in the 
increasing modulus of elasticity was determined. The reason for this may be the closure of 
micro and/or macro fractures under confining pressure. The mentioned finding is parallel 
to the findings obtained by Asef and Reddish. 

e In the rocks with high Schmidt hardness values, it was observed that the increase in the 
modulus of elasticity was higher when the confining stress was increased. However, the 
increase in the modulus of elasticity was found to be lower in the rocks with low strength 
and high porosity. 

e In the present research, a novel and practical approach to formulating the modulus of 
elasticity under three-axis loading conditions was applied. The suggested model is 
empirical, and experimental results constitute its basis. The verification of the GEP model's 
validity was performed by the statistical performance criteria utilized to assess the model, 
such as root mean square error (RMSE), mean absolute percentage error (MAPE), and the 
coefficient of determination (R’). These statistical parameters of GEP formulations also 
contributed to the excellent performance of GEP models in predicting the modulus of 
elasticity for natural stones. 

e The models proposed by Arora and Asef and Reddish exhibited a similar performance to 
each other in the current experimental data. However, the performance of these models is 
lower in comparison with the GEP-based model. 

e GEP is an appropriate method for estimating the EF value of rocks from anisotropic and 
heterogeneous materials as well as revealing nonlinear practical connections in which it is 
not possible to apply traditional techniques. Moreover, using GEP, the E value may be 
assessed without conducting complex and tiresome research facility tests. The proposed 
GEP formulations were easy since they could be used by anyone, who is not fully familiar 
with the GEP method, in comparison with other artificial intelligence methods. 
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